Warsztaty badawcze - praca domowa 2

LIME

Zbiór danych:

Wykorzystamy zbiór danych medycznych UCI Heart Disease, który zawiera wiek, płeć oraz wyniki badań medycznych pacjenta. Targetem jest ocena występowania wieńcowej choroby serca poprzez ocenę zwężenia naczyń wieńcowych (brak choroby - 0, choroba - 1). Zmienne kategoryczne (cp, thal oraz slope) zostały przetworzone za pomocą One-hot encoding, stąd w ramce danych pojawiły nam się zmienne z indeksami (np. thal_fd, thal_rd, thal_n).

Model:

Jako model wykorzystany zostanie Random Forest.

import pickle
import dalex as dx
import pandas as pd
import numpy as np
import warnings
warnings.filterwarnings('ignore')
from sklearn.model_selection import train_test_split

#wczytanie modelu
rf = pickle.load(open("./Modele/random_forest", 'rb'))

# wczytanie zbioru danych
data = pd.read_csv("./heart_data.csv")
data.head()

	age	sex	trestbps	chol	fbs	restecg	thalach	exang	oldpeak	ca	...	thal_fd	thal_rd	slope_up	slope_flat	slope_down	cp_ta	cp_aa	cp_np	cp_a	target
0	63	1	145	233	1	1	150	0	2.3	0	...	1	0	0	0	1	0	0	0	1	0
1	67	1	160	286	0	1	108	1	1.5	3	...	0	0	0	1	0	1	0	0	0	1
2	67	1	120	229	0	1	129	1	2.6	2	...	0	1	0	1	0	1	0	0	0	1
3	37	1	130	250	0	0	187	0	3.5	0	...	0	0	0	0	1	0	0	1	0	0
4	41	0	130	204	0	1	172	0	1.4	0	...	0	0	1	0	0	0	1	0	0	0

5 rows × 21 columns

data.info()

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 303 entries, 0 to 302
Data columns (total 21 columns):
 #   Column      Non-Null Count  Dtype  
---  ------      --------------  -----  
 0   age         303 non-null    int64  
 1   sex         303 non-null    int64  
 2   trestbps    303 non-null    int64  
 3   chol        303 non-null    int64  
 4   fbs         303 non-null    int64  
 5   restecg     303 non-null    int64  
 6   thalach     303 non-null    int64  
 7   exang       303 non-null    int64  
 8   oldpeak     303 non-null    float64
 9   ca          303 non-null    int64  
 10  thal_n      303 non-null    int64  
 11  thal_fd     303 non-null    int64  
 12  thal_rd     303 non-null    int64  
 13  slope_up    303 non-null    int64  
 14  slope_flat  303 non-null    int64  
 15  slope_down  303 non-null    int64  
 16  cp_ta       303 non-null    int64  
 17  cp_aa       303 non-null    int64  
 18  cp_np       303 non-null    int64  
 19  cp_a        303 non-null    int64  
 20  target      303 non-null    int64  
dtypes: float64(1), int64(20)
memory usage: 49.8 KB

Tak jak widzimy powyżej, nasz zbiór zawiera 21 kolumn. Zmienne thal, slope oraz cp zostały zakodowane za pomocą One-hot encoding.

# odzielenie targetu od innych zmiennych
y = data.target.values
x = data.drop(['target'], axis = 1)
x_train, x_test, y_train, y_test = train_test_split(x,y, test_size = 0.2,random_state=0, stratify=y)

# stworzenie explainera
explainer = dx.Explainer(rf, x_train, y_train)

Preparation of a new explainer is initiated

  -> data              : 242 rows 20 cols
  -> target variable   : 242 values
  -> model_class       : sklearn.ensemble._forest.RandomForestClassifier (default)
  -> label             : Not specified, model's class short name will be used. (default)
  -> predict function  : <function yhat_proba_default at 0x000002E5BAFD9D30> will be used (default)
  -> predict function  : Accepts pandas.DataFrame and numpy.ndarray.
  -> predicted values  : min = 0.0322, mean = 0.46, max = 0.988
  -> model type        : classification will be used (default)
  -> residual function : difference between y and yhat (default)
  -> residuals         : min = -0.738, mean = -0.00088, max = 0.81
  -> model_info        : package sklearn

A new explainer has been created!

Część 1: Predykcja i jej dekompozycja dla wybranej obserwacji.

Sprawdźmy jak działa nasz explainer w praktyce. Wybierzmy pierwszą obserwację w zbiorze danych oraz wyliczmy dla niej predykcję modelu.

x_train.iloc[0,:]

age            56.0
sex             1.0
trestbps      125.0
chol          249.0
fbs             1.0
restecg         1.0
thalach       144.0
exang           1.0
oldpeak         1.2
ca              1.0
thal_n          1.0
thal_fd         0.0
thal_rd         0.0
slope_up        0.0
slope_flat      1.0
slope_down      0.0
cp_ta           1.0
cp_aa           0.0
cp_np           0.0
cp_a            0.0
Name: 111, dtype: float64

Charakterystyka wybranego pacjenta (kilka wyróżniających się zmiennych):

• mężczyzna
• wiek 56 lat
• wysoki cholesterol (chol = 249.0)
• zmiany wskazujące na niedokrwienie mięśnia sercowego (restecg = 1)
• ból klatki typical angina często towarzyszący przy chorobie wieńcowej (cp_ta = 1)
• występowanie dławicy piersiowej (exang = 1)
• występowanie uwapnień (ca = 1)
• brak ubytków w badaniu scyntografii perfuzji mięśnia sercowego (thal_n = 1)

Po przyjrzeniu się danym możemy przypuszczać, że mężczyzna ten posiada chorobę wieńcową. Wartość targetu = 1 potwierdza nasze przypuszczenia.

y_train[0]

1

Predykcja modelu natomiast wynosi około 0.758.

explainer.predict(x_train)[0]

0.7583198920657599

Sprawdźmy dekompozycję predykcji za pomocą LIME:

ps = explainer.predict_surrogate(x_train.iloc[0,:], type = "lime", class_names=['zdrowy', 'chory'])
ps.show_in_notebook()

Wnioski:

• model zdecydowanie wskazuje na występowanie choroby, co zgadza się z wartością targetu.
• najważniejsze zmienne, zostały opisane już wyżej, ale warto podkreślić szczególnie istotność ca, cp, thal. Występowanie uwapnień i bólu klatki piersiowej typu typical angina wskazuje na chorobę, a brak ubytków w badaniu thal bardziej charakteryzuje osoby zdrowe.
• to co warto jeszcze zauważyć, a być może nie jest tak oczywiste, to wpływ zmiennej sex. Otóż okazuje się, że to u mężczyzn częściej występuje chorba wieńcowa, szególnie z objawami bólu klatki takiej jak u analizowanego pacjenta.

Na koniec jeszcze porównamy otrzymane z LIME wyniki z Shapley Values:

pp_shap1 = explainer.predict_parts(x_train.iloc[0,:], type='shap')
pp_shap1.plot()

Wnioski:

• wyniki podobne jak przy dekompozycji LIME - najważniejsze okazują się ca,cp,exang oraz thal.

Część 2: Dekompozycja LIME dla różnych obserwacji w zbiorze.

Zbadamy 5 obserwacji zarówno wskazujących na występowanie choroby, jak i jej brak. Sprawdzimy czy są jakieś zmienne, które mają duży wpływ na predykcje.

ps = explainer.predict_surrogate(x_train.iloc[13,:], type = "lime", class_names=['zdrowy', 'chory'], show_all=True)
ps.show_in_notebook()

ps = explainer.predict_surrogate(x_train.iloc[14,:], type = "lime", class_names=['zdrowy', 'chory'])
ps.show_in_notebook()

ps = explainer.predict_surrogate(x_train.iloc[15,:], type = "lime", class_names=['zdrowy', 'chory'])
ps.show_in_notebook()

ps = explainer.predict_surrogate(x_train.iloc[23,:], type = "lime", class_names=['zdrowy', 'chory'])
ps.show_in_notebook()

ps = explainer.predict_surrogate(x_test.iloc[54,:], type = "lime", class_names=['zdrowy', 'chory'])
ps.show_in_notebook()

Wnioski:

• zmienna ca - uwapnienie. Dla wartości 0 wskazuje osobę zdrową, a dla 1 (lub więcej) osobę chorą.
• thal_n - dla wartości 1 oznacza brak występowania ubytków i wskazuję na osobę zdrową. Gdy równy 0, oznacza, że występuje jedna z dwóch wersji ubtyków (stały - thal_fd, bądź odwracalny - thal_rd), które wskazują na chorobę.
• cp_ta - ból klatki piersiowej typical angina, który jeśli występuje to zdecydowanie wskazuje na chorobę.
• exang - dławica piersiowa. Jeśli występuje, to wskazuje na chorobę.
• sex - płeć. Gdy równa 1 (mężczyzna) możemy się spodziewać różnież dodatniej kontrubucji w predykcję choroby.
• oldpeak - obniżenie odcinka ST w elektrokardiograficznej próbie wysiłkowej. Dla dużych wartości (około > 1) wskazuje chorobę, dla mniejszych - jej brak.

Część 3: Podsumowanie.

• W każdej z analizowanych obserwacji wśród zmiennych mających największy wpływ na predykcje pojawiały się przede wszystkim cp,ca,thal. Są to zmiennie, które są dobrymi "wskaźnikami" występowania choroby wieńcowej serca. Posiadając dodatkowo wiedzę z zakresu medycyny i tego jakie te zmienne mają interpretacje, możemy poprzeć te stwierdzenie.
• wyjaśnienia LIME zgodne z naszą interpretacją danych oraz wyjaśnieniami Shapley Values